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Abstract 


We  investigate  the  linear  stability  of  static  coronal-loop  models 
undergoing  thermal  perturbations.  The  effect  of  conditions  at  the  loop 
base  on  the  stability  properties  of  the  models  is  considered  in  detail.  We 
consider  the  question  of  appropriate  boundary  conditions  at  the  loop  base 
and  conclude  that  the  most  physical  assumptions  are  that  the  temperature 
and  density  (or  pressure)  perturbations  vanish  there.  However,  if  the  base 
is  taken  to  be  sufficiently  deep  in  the  chromosphere,  either  several 
chromospheric  scale  heights  or  several  coronal  loop  lengths  in  depth,  then 
the  effect  of  the  boundary  conditions  on  loop  stability  becomes  negligible 
so  that  all  physically  acceptable  conditions  are  equally  appropriate.  For 
example,  one  could  as  well  assume  that  the  velocity  vanishes  at  the  base. 

We  calculate  the  growth  rates  and  eigenmodes  of  static  models  in  which 
gravity  is  neglected  and  in  which  the  coronal  heating  is  a  relatively 
simple  function,  either  constant  per-unit  mass  or  per-unit  volume.  We  find 
that  all  such  models  are  unstable  with  a  growth  rate  of  the  order  of  the 
coronal  cooling  time.  The  physical  implications  of  these  results  for  the 
solar  corona  and  transition  region  are  discussed. 


For  some  time  now  static  models  of  coronal  loops  (e.g. ,  Rosner, 
Tucker,  and  Vaiana  1978;  Craig,  McClymont,  and  Underwood  1978;  Vesecky, 
Actiochos,  and  Underwood  1979)  have  been  widely  used  to  interpret 
observations  of  both  solar  and  stellar  coronae  (e.g.,  Bonnett  and  Dupree 
1980,  Orrall  1981).  The  key  assumption  in  these  models  is  that  the  energy 
input  to  the  corona  is  constant  in  time  so  that  a  static  solution  to  the 
relevant  equations  is  possible.  However,  even  if  a  static  solution  is 
mathematically  possible,  it  will  not  be  physically  realizable  unless  this 
solution  is  also  stable  to  small  amplitude  perturbations.  It  is  well  known 
that  the  coronal  plasma  is  susceptible  to  a  thermal  instability  due  to  the 
form  of  the  dependence  of  radiative  loss  rate  on  temperature  (Field  1965). 
Therefore,  in  order  to  assess  their  degree  of  validity,  it  is  necessary  to 
determine  whether  the  static  models  are  thermally  stable  or  unstable. 

The  linear  theory  for  the  thermal  stability  of  static  coronal- loop 
models  has  been  investigated  by  a  number  of  authors,  but  with  differing 
conclusions.  In  his  original  work,  Antiochos  (1979)  concluded  that  the 
static  models  were  thermally  unstable.  Similar  results  were  obtained  by 
Hood  and  Priest  (1980).  However,  Chiuderi,  Einaudi,  and 
Torricelli-Ciamponi  (1981);  Craig  and  McClymont  (1981);  and  McClymont  and 
Craig  (1981a, b, c)  have  found  that  the  models  are  either  stable  or  that  the 
growth  rates  for  instability  are  too  small  to  be  physically  significant. 
The  origin  of  this  difference  in  the  results  of  the  two  sets  of  authors  is 
in  their  treatment  of  the  base  of  the  loop  models.  Antiochos  (1979)  and 
Hood  and  Priest  (1980)  have  not  included  cool  material,  T  K  10^  K,  in  their 
model,  for  which  the  form  of  radiative  loss  curve  (e.g.,  Raymond,  Cox,  and 
Smith  1976)  favors  linear  stability  (Field  1965).  In  addition,  Antiochos 


has  considered  only  perturbations  with  a  vanishing  first-order  heat  flux  at 
the  base.  Chiuderi,  Einaudi  and  Torricelli-Ciamponi;  and  Craig  and 
McClymont  argue  that  the  grcwth  rates  for  instability  are  very  sensitive  to 
these  assumptions  so  that  the  models  can  be  effectively  stabilized  either 
by  including  chromospheric  material  in  the  model  (Craig  and  McClymont  1981) 
or  by  changing  the  boundary  conditions  so  that  the  temperature 
perturbation,  instead  of  the  heat-flux  perturbation,  is  assumed  to  vanish 
at  the  base  (Chiuderi,  Einaudi,  and  Torricelli-Ciamponi  1981). 

The  reason  for  this  apparently  high  sensitivity  of  the  models  to  the 
base  conditions  lies  in  one  of  the  well-known  properties  of  the  static 
models:  the  fact  that  the  magnitude  of  the  conduction  and  radiation  terms 
in  the  energy  balance  are  approximately  equal  throughout  the  loop  (Vesecky, 
Antiochos,  and  Underwood  1978).  As  discussed  by  Field  (1965),  we  expect 
thermal  instability  in  coronal  plasma  whose  structure  is  such  that 
radiation  dominates.  On  the  other  hand,  we  expect  stability  if  the 
structure  is  conduction  dominated,  since  conduction  always  acts  to  damp  out 
any  temperature  perturbation.  Hence,  the  static  models,  whose  structure  is 
such  that  these  two  terms  are  comparable,  are  perched  near  the  boundary 
separating  stability  from  instability.  Even  seemingly  minor  changes  on  the 
form  of  the  base  or  of  the  boundary  conditions  at  the  base  can  be 
sufficient  to  push  the  models  to  one  side  or  the  other  of  this  boundary. 

It  is  clear,  therefore,  that  the  proper  treatment  of  the  base  is 
critical  for  determining  to  what  extent  the  models  are  either  stable  or 
unstable.  In  the  next  section  we  discuss  the  question  of  the  base 
conditions  in  detail.  In  Sections  III  and  IV  we  calculate  growth  rates  and 
eigenmodes  for  various  static  models  and  for  various  boundary  conditions. 
In  the  final  section,  V,  we  discuss  the  implications  of  our  results. 


II.  Boundary  Confll tiflflg 


As  an  aid  to  determining  the  appropriate  boundary  conditions,  let  us 
first  write  down  the  perturbation  equations  for  a  one-dimensional  loop 
model.  Assuming  that  all  physical  variables  are  of  the  form: 

f(s,t)  =  fQ(s)  +  eVt  f^s)  ,  (1) 
and  assuming  that  the  growth  or  damping  time  1/v  is  slew  compared  to  the 
sound  travel  time  so  that  the  acceleration  terms  can  be  neglected  in  the 
momentum  equation,  we  obtain  (Field  1965,  Antiochos  1979): 


vn.  +  -7  (An  V, )  =  0 
1  A  ds  0  1 


(2) 


d_ 

ds 


m  n1  g 


(3) 


3/2<«Pl  +  V  P0)  -  5/2  *0  (vn  +  V  fj  n  ) 

"0 


+  i  s  ■  £i  w 


and 


(5) 


where:  A(s)  is  the  cross-sectional  area  of  the  loop;  g|j  (s)  is  the 
component  of  gravity  parallel  to  the  loop;  is  the  perturbed  heat  flux, 
which  for  the  Spitzer  (1962)  conductivity  is  given  explicitly  by: 


10-6  —  (T  T  ) 
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constitutes  the  perturbed  energy  sources  and  sinks  of  the  plasma.  For 
optically  thin  radiative  losses  and  for  a  coronal  heating  function  that 
depends  only  on  density  and  temperature  or,  equivalently,  pressure  and 
temperature,  is  given  explicitly  by  (Antiochos  1979): 
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where  A (T)  is  the  radiative  loss  coefficient  (e.g. ,  Cox  and  Tucker 
and  E(p,T)  is  the  coronal  heating  rate.  The  equilibrium  profiles 
course,  given  by  (e.g.,  Vesecky,  Antiochos,  and  Underwood  1979): 


(7) 
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Note  that  we  have  neglected  in  equations  (2)  -  (8)  effects  due  to 
changes  in  ionization  equilibrium  or  to  optically  thick  radiation,  both  of 
which  are  clearly  important  in  the  chromosphere.  Hence,  this  model  cannot 
be  used  to  determine  the  stability  of  the  true  solar  chromosphere.  We  use 
this  model  chromosphere  only  for  investigating  the  effect  of  having 
thermally  stable  material  at  the  base  of  the  corona  and  transition  region. 
Since  our  results  will  turn  out  to  be  insensitive  to  the  detailed  structure 
of  this  region,  we  expect  that  our  analysis  will  also  be  applicable  to  the 
solar  chromosphere  as  long  as  it  is  thermally  stable  as  well. 

It  is  evident  from  the  equations  above  that  there  are  four  independent 
spatial  derivatives  in  the  problem,  which  implies  that  four  boundary 
conditions  are  required.  This  is  to  be  expected  since  the  full  nonlinear 
equations  also  require  four  spatial  boundary  conditions  for  a  unique 
solution  (e.g.,  Richtmyer  and  Morton  1967).  The  usual  situation  is  that  at 
end  of  the  loop  two  boundary  conditions  are  specified.  One  condition 
generally  relates  to  the  thermal  properties  of  the  base;  for  example,  with 
respect  to  the  perturbation  the  loop  base  may  act  as  a  thermal  bath,  in 
which  case  the  appropriate  condition  is  that  the  temperature  perturbation 
vanishes  there,  Tlb  =  0;  or  it  may  act  as  a  thermal  insulator  so  that  the 


perturbation  heat  flux  vanishes  there,  F,,  =  0.  The  other  condition 

lb 

generally  relates  to  the  inertial  properties  of  the  base;  for  example,  it 
may  act  as  a  rigid  wall  with  respect  to  the  perturbation  so  that  the 
velocity  vanishes,  Vlb  =  0,  or  a  "free  surface,"  so  that  the  pressure 
perturbation,  P  =  0.  Most  of  the  discussion  in  the  literature  has 
concerned  these  four  conditions,  and  the  discussion  in  this  paper  will  also 
concentrate  on  these  four.  Let  us  emphasize,  however,  that  there  are  an 
infinite  number  of  other  possible  boundary  conditions  and  that,  lacking 
some  physical  justification,  there  is  no  reason  to  single  out  the  four 
above.  A  priori,  all  conditions  are  equally  acceptable. 

There  has  also  been  a  considerable  discussion  in  the  literature  on  the 
distinction  between  symmetric  and  antisymmetric  modes.  However,  this 
distinction  is  spurious.  It  is  due  only  to  the  special  class  of  models 
considered  by  the  previous  authors  (Antiochos  1979;  Hood  and  Priest  1980; 
Chiuderi,  Einaudi,  and  Torricelli-Ciamponi  1981;  Craig  and  McClymont  1981). 
In  their  models  gravity  is  neglected,  and  the  loop  area  and  the  coronal 
heating  do  not  have  any  spatial  variation.  Under  these  simplifications, 
the  static  equations  are  autonomous,  and  the  equilibrium  models  are 
symmetric  about  the  loop  apex.  Hence,  the  solutions  to  the  first-order 
equations,  i. e. ,  the  normal  modes  of  the  loop,  are  either  purely  symmetric 
or  antisymmetric  about  the  apex.  These  results  will  clearly  not  hold  for  a 
realistic  equilibrium  model  that  has  no  special  symmetry  properties.  For 
such  a  model  each  choice  of  four  boundary  conditions  will,  in  general, 
determine  a  complete  and  distinct  set  of  eigenvalues  and  eigenfunctions. 
Within  each  set  all  the  modes  are  equally  valid.  Hence,  in  this  paper  we 
will  avoid  any  discussion  in  terms  of  symmetric  or  antisymmetric  modes;  and 
since  the  most  important  mode  from  the  viewpoint  of  thermal  instability  is 
the  lowest  order  one,  we  will  concentrate  on  the  fundamental. 


The  question  remains  as  to  what  choice  of  boundary  conditions  best 
represents  the  physical  situation  at  the  base  of  a  solar  loop.  Ideally, 
one  would  like  conditions  at  the  loop  base  not  to  affect  the  stability 
properties  of  the  corona  and  transition  region.  The  standard  procedure  for 
achieving  this  is  to  place  the  base  sufficiently  far  that  over  the  time 
scale  of  any  instability  no  effects  can  propagate  to  the  base.  However, 
the  sound  travel  time  in  the  corona,  chromosphere  and  even  photosphere  is 
very  rapid  compared  to  typical  thermal  time  scales,  such  as  the  coronal 
cooling  time.  On  this  time  scale,  material  at  the  loop  base  is  always  in 
contact  with  the  coronal  plasma,  even  if  the  base  is  assumed  to  be  in  the 
photosphere.  Note  that  our  neglect  of  the  acceleration  terms  in  the  force 
equation  is  tantamount  to  assuming  an  infinite  sound  speed  and, 
consequently,  a  negligible  propagation  time  scale  throughout  the  loop.  For 
such  a  model  it  is  not  obvious,  at  least  to  us,  whether  the  position  of  the 
base  can  ever  be  chosen  so  that  conditions  there  have  an  insignificant 
effect  on  the  behavior  of  the  corona.  We  show  below  that,  fortunately, 
even  with  the  rapid  sound  speed,  the  base  has  a  negligible  influence  on  the 
coronal  stability  if  it  is  placed  sufficiently  deep. 

In  order  to  determine  the  physically  valid  boundary  conditions  to 
impose  at  the  loop  base,  we  need  a  definite  model  for  this  base  region.  As 
a  simple  model  let  us  assume  that  the  chromosphere  below  a  coronal  loop  is 
isothermal  at  temperature  Tb  ,  is  in  hydrostatic  equilibrium,  and  is 
strictly  stable.  In  addition,  we  assume  there  is  no  real  end  point  to  a 
solar  loop,  so  this  model  base  region  extends  below  the  corona  and 
transition  region  to  arbitrarily  large  depths,  i.e. ,  either  to  many 
chromospheric  scale  heights  or  to  depths  large  compared  to  the  coronal  loop 
length.  Finally,  let  us  neglect  the  curvature  of  the  loop  for  this 
analysis  and  set  g u  (s)  to  be  a  constant. 


The  form  of  the  boundary  conditions  follows  simply  from  the 
requirement  that,  since  we  are  interested  in  coronal  disturbances,  the 
perturbations  should  vanish  or  at  least  remain  finite  in  the  chromosphere. 
Combining  equations  (3)  and  (5),  an  expression  can  be  obtained  for  the 

pressure  perturbation  P  in  terms  of  the  temperature  perturbation  T  : 

s 

fT_ 1A.  (9) 

p0(s)  P0<0)  J  T0  H(T0)  ’ 

where  H  (Tq)  is  the  gravitational  scale  height  at  temperature  T^: 

H (T)  s  — (10) 

m8H 

and  s  =  0  is  some  reference  point  in  the  loop  (the  top  is  convenient).  We 

are  defining  s  to  increase  downwards.  It  is  evident  from  (9)  that  in  order 

for  P^s)  to  remain  finite  with  increasing  depth  in  the  chromosphere,  i.e. , 

s  ->  °°,  the  temperature  perturbation  must  vanish  with  depth  ->  0.  In 

fact,  must  decay  quite  rapidly  with  depth,  on  the  size  scale  of  the 

chromospheric  scale  height  H(T,  )  or  less,  because  the  equilibrium  pressure 

b 

PQ  grows  exponentially  on  this  size  scale.  Note  that  we  are  excluding  the 
possibility  of  modes  that  oscillate  indefinitely  with  depth.  We  are 
interested  primarily  in  the  lowest  order  mode;  and  as  we  will  see  below, 
this  usually  has  no  zero  crossings. 

We  conclude,  therefore,  that  the  proper  thermal  boundary  condition  is 

that  T  =  0.  This  result  is  to  be  expected  physically.  It  follows  from 

lb 

the  fact  that  any  perturbation  of  the  temperature  also  implies  a 
perturbation  of  the  gravitational  scale  height.  If  the  temperature 
disturbance  extends  more  than  a  scale  height  or  so  deep  into  the 
atmosphere,  a  very  large  pressure  disturbance  must  result  due  to  the  change 
in  scale  height.  Hence,  the  only  physically  acceptable  temperature 
perturbations  are  those  that  vanish  on  the  order  of  a  scale  height  or  less. 


Now,  to  determine  the  inertial  boundary  condition,  we  integrate  the 


continuity  equation  (2)  to  obtain  an  expression  for  the  particle  flux  n 
in  terms  of  the  density  perturbation  n^: 

s 

n  (s)V  (s)  =  n  (0)  V. (0)  -  vl  n  ds  .  (11) 

0  1  0  1  11 
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As  in  the  situation  above,  it  is  evident  that  in  order  for  to  remain 

finite  as  s  ->  00 ,  n^  must  vanish  with  depth.  Again,  this  is  to  be  expected 
physically.  If  the  density  disturbance  extends  over  a  large  depth,  then 
the  change  in  the  total  number  of  particles  in  this  region  must  be  large, 
and  consequently  a  large  flux  of  particles  into  or  out  of  the  region  must 
result.  It  is  evident  from  (11)  that  in  the  chromosphere  n^  must  vanish 
faster  than  Lc/s,  where  Lcis  the  size  scale  of  the  coronal  portion  of  the 
loop,  i.e. ,  the  region  where  we  expect  a  finite  n..  Note  that  the  density 
perturbation  may  extend  to  much  greater  depths  than  the  temperature  one. 
We  found  above  that  must  vanish  at  least  as  fast  as  expf-s/I^),  where  1^ 
is  the  gravitational  scale  height  in  the  chromosphere.  For  typical  solar 
loops  1^  «  lc. 

We  conclude  that  the  proper  boundary  conditions  at  the  loop  base  are 
that  Tlb  =  0  and  nlb  =0.  At  any  finite  depth  this  is  equivalent  to  the 
set  Tlb  =  0  and  Plb  =  0.  Although  the  temperature,  density,  and  pressure 
perturbations  vanish,  the  mass  flux  does  not.  Using  (2)  and  (3),  the 
perturbed  mass  flux  can  be  related  to  the  pressure  perturbation: 


novi  = 


n0(°) 


V1(0) 


(PjGO  -  p  (0)) 


(12) 


Since  the  pressure  perturbation  vanishes  with  depth,  the  mass  flux  will  in 
general  tend  to  some  finite  constant,  corresponding  to  a  steady-state  flow 
in  the  chromosphere  and  below.  This  is  to  be  expected  on  physical  grounds 
as  well.  We  do  not  expect  the  temperature  and  density  structure  deep  in 


the  chromosphere  or  photosphere  to  be  affected  by  a  disturbance  in  the 
corona.  However,  if  as  in  our  problem,  the  time  scale  of  the  disturbance 
is  long  compared  to  the  sound  travel  time,  then  a  steady  flow  must  occur  at 
large  depths  to  support  any  coronal  flows.  Of  course,  the  velocity  itself, 
Vj  ,  vanishes  with  depth  since  the  density  increases;  however,  the  mass 
flux  stays  constant.  This  is  exactly  the  situation  that  we  found  for 

the  case  of  evaporative  cooling  of  coronal  loops  (Antiochos  and  Sturrock 
1978).  There  again,  the  coronal  evolution  occurred  on  time  scales  long 
compared  to  the  sound  travel  time,  and  it  was  found  that  a  steady  flow  was 
set  up  at  the  loop  base  due  to  evaporative  motions  in  the  corona.  This  is 
also  the  situation  for  the  photosphere  at  the  base  of  an  open  magnetic 
field  region  such  as  in  a  coronal  hole.  In  these  regions  coronal  heating 
results  in  a  steady  solar  wind  flew.  The  presence  of  this  flow  has  no 
observable  effect  on  photospheric  temperatures  and  pressures,  but  it 
clearly  implies  the  existence  of  a  mass  flux  in  the  photosphere  equal  to 
the  coronal  mass  flux. 

Let  us  now  examine  the  effect  of  using  different  boundary  conditions 
at  the  loop  base,  in  particular,  the  rigid-wall  one.  We  wish  to  find  out 
how  the  coronal  perturbation  depends  on  the  base  conditions,  for  example, 
how  P,  (0)  depends  on  n„.  V1t  .  From  (12)  it  is  evident  that  specifying  the 
base  mass  flux  is  equivalent  to  specifying  the  base  pressure  perturbation, 


but  the  relation  between  P. (0)  and  P_,  is  given  by  (9): 

1  lb  L 


V°>  =  P 


P!b  +  P0(0) 


11  ds 
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Equation  (13)  is  the  key  result.  The  effect  on  the  coronal  perturbation 


Pj  (0)  of  assuming  a  non- vanishing  P^  is  given  by  the  first  term  on  the 
right-hand  side  of  (13).  This  term  is  of  order  P  (0)/P  for  a  temperature 


perturbation  that  vanishes  rapidly  in  the  chromosphere  so  that  the  pressure 


C  =10  and  m  =  3.  This  value  for  m  is  convenient  numerically  since  for 


m  =  3  the  depth  of  the  base  region  turns  out  to  be  proportional  to 
1/6;  hence,  the  range  of  chromospheric  depths  covered  by  the  range  used  for 
£  is  correspondingly  large,  (L  -  Lc)/Lc  ranges  from  “  10  5  to  10.  Again, 
we  note  that  there  is  almost  no  variation  in  v  irrespective  of  the  depth 
of  the  chromosphere. 

We  have  investigated  the  effect  of  the  form  of  the  coronal  heating  on 
stability  by  performing  a  series  of  calculations  as  in  Figure  3  for  the 
case  of  uniform  heating  per  unit  volume,  Y  =  C.  The  results  are 
essentially  identical  to  those  shown  in  Figure  3.  This  agrees  with  our 
previous  calculations  (Antiochos  1979)  in  which  we  investigated  a  wide 
range  of  values  for  Y  and  found  no  significant  effect  on  the  growth  rate. 

We  have  also  looked  for  any  possible  effect  due  to  the  form  of  the 
chromospheric  heating  by  varying  the  exponent  m  defined  in  (43).  This 
parameter  determines  the  degree  of  stability  of  the  chromosphere,  i.e. ,  the 
magnitude  of  the  quantity  ( T0'--f/3  TQ  )  (cf . ,  Field  1965).  The  chromospheric 
magnitude  of  this  quantity  is  small  in  our  models  compared  to  its  coronal 
value  except  for  the  case  m  =  1.  Note  that  for  m  =  1  the  chromospheric 
depth  L  -  L  varies  only  as  ln(1/6  )  so  that  even  with  the  smallest  value  of 
6  usable  numerically,  ~  10  ^ ,  the  chromospheric  depth  is  not  large, 
(L  -  Lc)/Lc  ”  .05.  The  results  for  the  cases  with  m  =  1  are  identical  to 
those  shown  in  Figure  3  for  the  case  m  =  3.  Again,  there  is  no  significant 
dependence  of  the  growth  rate  v  on  either  the  base  temperature  or  the 
chromospheric  depth.  We  conclude  that  the  instability  does  not  depend  on 

the  particular  form,  (42),  that  we  assumed  for  the  chromospheric  heating. 

4 

This  result  is  to  be  expected.  For  temperatures  K  5  x  10  K,  the 
quantity  T is  strictly  negative  in  all  our  models;  hence,  this  region 
cannot  contribute  to  any  instability.  As  we  will  see  belcw,  the  unstable 


obtained  by  solving  (14)  for  a  series  of  static  models  with  different  base 
properties.  These  static  models  all  have  very  similar  coronal  properties. 

The  apex  temperature  T^(0)  was  chosen  in  each  case  to  be  106  K  and  the 

-12 

coronal  heating  given  by  y  =  1  and  Q  =  10  ergs/see/partiele,  equation 
(41).  Given  the  coronal  temperature  and  heating  rate,  then  as  is  well 
known  (e.g. ,  Vesecky,  Antiochos,  and  Underwood  1979),  the  coronal  density 
and  coronal  loop  length  are  constrained  by  the  scaling  laws.  They  turn  out 

q  —  3  8 

to  be  n„(0)  =  3.2  x  10  cm  and  L  =  5.2  x  10  cm.  Conditions  at  the 

u  c 

base,  however,  are  very  insensitive  to  the  coronal  conditions,  so  that  we 
can  vary  the  base  temperature  and  the  amount  of  base  material  with  almost 
no  change  in  coronal  parameters. 

In  our  model  the  base  temperature  is  given  by  Tb,  and  the  depth  of  the 
chromosphere  is  determined  primarily  by  the  constant  6,  equation  (43). 
Hence,  there  are  actually  two  distinct  sets  of  models  in  Figure  3.  In  one 
set  all  the  models  have  a  value  of  6=1,  which  implies  that  there  is  an 
insignificant  amount  of  material  at  the  base.  These  models  correspond  to 
the  usual  static-loop  models,  as  in  Vesecky,  Antiochos,  and  Underwood 

(1979).  For  such  models  we  investigated  the  effect  of  the  base  temperature 

5  4 

on  stability  by  varying  T(j  from  2  x  10  K  to  10  K.  In  Figure  3  only  the 

4 

cases  down  to  T^  =  3  x  10  K  are  shown,  but  no  change  was  observed  in  using 

4 

T  down  to  10  .  As  is  clear  from  the  figure,  the  value  of  the  base 
temperature  has  essentially  no  effect  on  loop  stability.  The  growth  rate 
vr  is  of  order  unity  irrespective  of  whether  T^  is  chosen  to  occur  at 
chromospheric  or  transition-region  temperatures. 

The  second  set  of  models  in  Figure  3  are  those  with  a  fixed  base 

4 

temperature  Tfa  =  3  x  10  K,  but  with  varying  amounts  of  base  material  at 

— r> 

this  temperature.  The  value  of  £  was  varied  from  ?  =  1  to  '  =5x10 

The  other  parameters  in  the  energy  input  rate  (equation  (42))  were  fixed, 


such  as  the  VAL  models  (Vernazza,  Avrett  and  Loesser  1981),  all  indicate 

4 

the  presence  of  a  flat  temperature  "plateau"  at  T  "  2  x  10  K  with  a  very 
abrupt  rise  at  higher  temperatures.  Within  the  context  of  the  static 
model,  and  neglecting  possible  kinetic  effects  (Shoub  1983),  the  only  way 
to  obtain  the  observed  very  sharp  structure  in  the  temperature  profile  is 
to  have  a  correspondingly  sharp  structure  in  the  energy  input  rate.  This 
is  one  reason  why  the  static  models  are  generally  assumed  to  have  a  base 
above  2  x  10^K.  The  particular  form,  (41)  and  (42)  is  chosen  simply  for 
convenience.  It  permits  us  to  investigate  the  effects  of  adding  a  varying 
amount  of  cool,  radiatively  stable  material  to  the  loop  base. 

b)  Growth  Rates 
(i)  Free  Surface 

We  have  calculated  the  eigenvalues  A  for  a  wide  range  of  static 
models.  In  particular,  we  have  determined  the  dependence  of  the  stability 
on  the  following  properties  of  the  model:  the  coronal  temperature,  the 

base  temperature,  the  depth  of  the  chromosphere,  the  form  for  the  heating 
and,  of  course,  the  boundary  conditions.  The  main  result  of  this  analysis 
is  that  for  the  boundary  conditions  T^b  =  0,  n^b  =  0,  all  the  models  were 
unstable  with  a  growth  time  approximately  equal  to  the  coronal  cooling 
time.  The  growth  rate  was  found  to  be  completely  Insensitive  to  the  loop 
temperature  either  in  the  corona  or  base,  to  the  amount  of  chromospheric 
material  or  to  the  form  of  the  chromospheric  heating.  Hence,  we  conclude 
that  the  instability  is  a  true  physical  one  and  will  occur  in  solar  loops 
that  satisfy  our  basic  assumptions,  i.e. ,  spatial  variations  in  the  loop 
area  and  heating,  and  the  effects  of  gravity  can  be  neglected. 

Figure  3  shews  the  growth  rate  vt  of  the  lowest  free-surface  mode 


T  <  5  x  104  K.  In  order  for  the  model  to  have  a  large  chromosphere,  the 
temperature  gradients  and,  hence,  the  heat  flux  must  be  small  at  these 
temperatures.  The  radiation  losses  in  the  chromosphere  cannot  be  balanced 
by  conduction  from  above  (e.g. ,  Athay  1981);  thus  for  a  static  model  these 
losses  must  be  balanced  by  the  energy  input  (e.g.,  Craig,  Robb  and  Rollo 
1982). 

Therefore,  we  use  the  following  form  for  the  energy  input: 

c  =  Qn^  tanh(x  (T) )  +  n^A(1  -  tanh(x(T)))  ,  (41) 


where 


X (T)  =  C  (T/Tb  -  1  +  6)m 


(42) 


and  Q,y,  C  ,  m  and  6  are  constants.  Equations  (41)  and  (42)  imply  that  for 

6«  1: 


e  -> 


for  T  >  T, 

D 

for  T  =  T. 

b 


(43) 


V 

i.e. ,  the  energy  input  rate  tends  to  Qn  for  temperatures  larger  than  the 
base  temperatures  and  tends  to  the  radiative  loss  rate  at  the  base 
temperature.  The  values  0  and  1  for  y  correspond  to  constant  coronal 
heating  per  unit  volume  and  per  unit  mass,  respectively.  The  depth  of 
material  at  the  base  temperature  is  determined  by  the  value  of  5;  it 
becomes  arbitrarily  large  for  arbitrarily  small  6. 

The  form  that  we  assume  for  the  energy  input,  (41)  and  (42)  may  appear 
highly  unphysical  because  it  undegoes  an  abrupt  although  continuous  change 
near  T  =  T^.  However,  such  a  sharp  structure  in  the  energy  input  rate  is 


implied  by  the  observations.  Empirical  models  of  the  upper  chromosphere, 


chromosphere.  This  makes  it  impossible  to  use  a  uniform  grid  for  finite 
differencing  the  independent  variable  x.  However,  we  cannot  use  the 
equilibrium  temperature  as  the  independent  variable  since  the 

temperature  gradient  vanishes  at  the  loop  apex  and  base.  Instead,  we 
define  a  new  independent  variable  r: 


and  use  a  uniform  grid  in  r.  It  is  evident  from  (30)  that  the  effect  of 
transforming  to  r  is  to  increase  the  number  of  grid  points  at  transition 
region  temperatures  (TQ(0)  <  TQ<  Tfa)  by  a  factor  of  order  (TQ  (0)/TQ)2. 
We  found  that  this  was  sufficient  to  yield  an  accurate  solution. 

a)  Form  aL  Heating  and  Coaling 

There  are  two  functions  in  the  problem  that  have  yet  to  be  specified: 
the  radiative  loss  coefficient,  A(T),  and  the  energy  input  rate  c(n, T). 
For  the  former  we  use  a  smooth  analytic  approximation  to  the  curve  derived 
by  Raymond,  Cox  and  Smith  (1976).  Figure  2  shows  this  curve  along  with  our 
approximation  to  it.  The  stability  properties  of  the  models  are  not 
sensitive  to  the  exact  form  of  A(T)  as  long  as  this  form  is  such  that  it 
implies  radiative  instability  above  ”  10 5  K  and  stability  belcw. 

Since  the  mechanism  for  coronal  heating  is  not  known,  the  form  of  the 
energy  input  e  is  essentially  arbitrary.  In  our  previous  work  we  found 
that  the  stability  properties  are  insensitive  to  the  exact  form  for  e,  at 
least,  for  the  case  where  this  form  is  a  power-law  dependence  on  n  and  T. 
Hence,  we  simply  assume  that  in  the  corona  and  transition  region  the  energy 
input  either  per  unit  volume  or  per  particle  is  constant.  At  lcwer 
temperatures  a  different  form  for  the  energy  input  must  be  used  in  order  to 
have  a  significant  mass  of  material  at  "chromospheric"  temperature, 


models  the  growth  time  1  A>f  is  very  short,  of  order  the  cooling  time  in  the 
transition  region  rather  than  that  in  the  corona  (Antiochos  1979).  We  find 
that  for  these  cases  the  rigid-wall  growth  rate  vr  is  slightly  larger  than 
Vj.  However,  for  models  in  which  some  cool  stable  material  occurs  at  the 
base,  the  rigid-wall  boundary  conditions  always  result  in  smaller  growth 
rates. 

IV.  Numerical  Results 

The  mathematical  problem  that  must  be  solved  is  now  clearly  defined. 
For  the  free-surface  conditions,  (14)  reduces  to  a  standard  Sturm-Liouville 
problem,  as  in  Antiochos  (1979),  which  we  solve  in  exactly  the  same  manner 
as  before. 

For  the  rigid-wall  conditions,  the  inhomogeneous  problem  must  be 
solved  subject  to  the  constraint  (24a).  Note  that  this  is  not  a  standard 
Sturm-Liouville  problem;  hence,  none  of  the  well-known  theorems  (e.g. , 
Morse  and  Feshbach  1953)  on  the  behavior  of  the  eigenfunctions  and 
eigenvalues  need  apply.  There  are  two  parameters,  A  and  n ,  that  must  be 
determined;  therefore,  we  do  not  use  a  "shooting"  technique  as  in  the 
homogeneous  case  to  obtain  solutions.  First,  we  note  that  since  n  will,  in 
general,  be  non-zero  for  any  finite- length  loop,  we  can  eliminate  it  by 
defining  a  new  dependent  variable  (y/n).  Now  (14)  and  (24a)  define  a 
nonlinear  boundary- value  problem  for  (y/y)  and  the  unknown  parameter  A. 
This  is  solved  by  standard  techniques:  Newton-Raphson  iteration  and 
finite-difference  solution  of  the  linearized  equations.  Convergence  was 
rapid,  <  20  iterations  for  all  the  cases  we  investigated. 

One  numerical  complication  is  the  extreme  variation  in  temperature 
scale  height  between  the  corona,  the  transition  region,  and  the 


order  mode,  has  no  zero  crossings  and  its  normalization  is  irrelevant,  the 
sign  of  6v  is  determined  primarily  by  the  sign  of  the  terms  in  the  brackets 
in  the  rightmost  integral;  however,  this  factor  does  not  have  a  unique 
sign.  The  first  term  v  is  always  positive  and  for  models  of  interest  has 
magnitude  ”  1/x;  consequently,  it  enhances  instability.  The  sum  of  the 
second  and  third  term  is  proportional  to  (-e/P^ )  for  the  case  where  the 
heating  c  is  constant  per  particle,  and  (-2e/P^  for  e  constant  per  unit 
volume.  Hence,  the  contribution  from  the  last  two  terms  is  negative,  and 
since  in  the  static  model  the  energy  input  and  the  radiative  losses  are 
roughly  equal  in  the  corona  (e.g. ,  Vesecky,  Antiochos,  and  Underwood  1979), 
e/Pq  is  also  of  order  1/r.  It  is,  therefore,  not  clear  which  of  the  two 
contributions  will  dominate;  however,  it  turns  out  that  the  negative 
contribution  is  usually  larger  than  the  positive  one  so  that  the  rigid  wall 
has  a  net  stabilizing  effect. 

This  is  to  be  expected  physically.  A  rigid  wall  acts  to  inhibit  flows 
and,  consequently,  any  density  disturbance.  For  optically  thin  radiation 
at  coronal  temperatures,  both  the  temperature  dependence  and  the  density 
dependence  of  the  losses  is  such  as  to  promote  instability  (Field  1 96 5 ) . 
By  suppressing  the  density  effect,  a  rigid  wall  tends  to  weaken  the 
instability.  It  is  interesting  to  note  that  this  need  not  always  be  the 
case.  Sometimes  a  rigid  wall  can  actually  enhance  instability.  In 
particular,  if  v  is  very  large  and  dominates  the  other  terms  in  (39),  we 
would  expect  6v  to  be  positive.  This  situation  does  occur  for  models  in 
which  the  base  is  placed  high  up  in  the  transition  region  and  the 

temperature  boundary  condition  is  assumed  to  be  dT  /ds  =  0.  For  such 

lb 


the  particular  solution: 


2  5/14  (Vr  Tb  ~  5/2) 

5  nr  yb  (v  t  +  1  ) 

r  b 


If,  as  expected,  v 


“  0(  1 A  ) ,  then  for  tA  ->  0 
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2  5/14 

^  5  \  yb 


We  will  see  below  that  (35)  also  agrees  with  exact  solution  of  (14). 

These  results  verify  our  claim  that  n  decays  as  L  /L.  We  have  that 


c 

hr 
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-5/14  j  . 
y  dx  + 


h' 


-5/14  . 
y  dx 


Substitutingithe  asymptotic  form  (35),  or  (36) ,  into  the  second  integral 
yields  (26)  to  within  a  factor  of  order  unity. 

The  effect  on  the  growth  rate  of  imposing  a  rigid-wall  constraint  can 
also  be  estimated.  Using  the  equations  for  5  and  C  ,  and  their  boundary 

conditions,  the  following  expression  is  obtained: 

1  1. 


6A  =  A  -  A 


r  f  1 


f  f  f(x)  dx  +  f  £f  g(x)  dx 

i) _ 2 _ 

/  ^ d* 


If  L£/L  «  1,  then  to  lowest  order  in  Lc/L,  Equation  (38)  can  be  written 


6v  =  vr  -  vf  =1  — 


r  -5/14  \  1 

Cf  y  dxl  ( 


2 

dx 


!(vf-  2Po  +  i^)dx  •  (39) 


Each  of  the  integrals  in  (39)  is  of  order  Lc/L;  hence,  6v  is  also  of  this 


order. 


Depending  on  the  sign  of  6v,  the  rigid-wall  can  have  a  stabilizing  (6v 


negative)  or  destabilizing  (<$v  positive)  influence.  Since  t, f,  the  lowest 


The  longest  length  scale  for  the  decay  is  obtained  in  the  limit  t/t,  ->  0. 

b 

7  /A 

This  length,  "  Lc(Tb/TQ(0))  ,  is  typically  small  compared  to  physical 
size  scales  either  in  the  corona  or  chromosphere.  We  show  below  that  the 
form  predicted  above  is  obtained  from  exact  solution  of  (14). 

We  now  consider  the  rigid-wall  case.  With  £  as  the  dependent 
variable,  Equation  (14)  becomes: 


d_ 

dx 


(y  Cr)  +  Q(x)  +  Ar  =  \  (*r  f(x)  +  g(x))  1 


(30) 


where  the  subscript  "r"  is  used  to  indicate  rigid-wall  quantities.  The 
functions  f  and  g  are  given  by: 


ft  v  2  5/14 

f(x)  =  y  y 


(3D 


and 


g(x)  =-y 


5/14 


+ 

dx2 


_  ^ 

7.2  T«p0^ 


(32) 


Following  the  arguments  above,  in  the  base  region,  Equation  (30)  reduces 
to:  .2 


f-  r.  -  (k-Y  ,  v  r 

'  \\)  ’l.  r  \Lc/  r  r 

\  (r)  (-  !  \  T  *  T/Tb)  • 


(33) 


where,  again,  due  to  the  form  of  the  radiative  losses,  we  assume  that. 


( 


3X 

3P 


(34) 
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The  general  solution  to  (33)  consists  of  a  linear  combination  of  the 
homogeneous  and  particular  solutions.  However,  the  homogeneous  solution  is 
equivalent  to  the  free-surface  solution,  and  we  found  above  that  this 
decays  exponentially.  Hence,  the  solution  at  large  depths  is  dominated  by 


where  the  brackets  indicate  the  average  coronal  value.  Hence,  for  large 
depths,  n  decays  as  Lc/L.  This  is  the  same  behavior  that  we  found  before, 
except  that  now  since  there  is  no  gravity,  the  only  scale  in  the  problem  is 
the  coronal  loop  length,  so  that  this  length  rather  than  the  gravitational 
scale  height  determines  the  decay  scale  for  the  effects  of  the  wall. 


b)  Asymptotic 

In  order  to  check  whether  C  will,  in  fact,  decrease  rapidly  in  the 
chromosphere,  let  us  examine  the  asymptotic  forms  for  t,  in  this  region  as 
predicted  by  equation  (I1*).  We  first  consider  the  free-surface  case. 
Since  is  approximately  constant  in  the  chromosphere.  Equation  (14) 


reduces  to: 


dx 


2  ?f  + 


"(To  sfj  -  (lJ 


vfT5f 


(27) 


where  the  subscript  "f"  is  used  to  refer  to  the  free-surface  variables. 
Assuming  that  the  radiative  loss  function  has  a  simple  power-law  dependence 
on  T,  and  that  the  power  is  positive  so  that  the  chromosphere  is 
radiatively  stable,  then  we  expect  that: 


lT_o  ax  \  i_ 
\po9To /  '  Po 


-  1/t, 


(28) 


where  t.  is  the  radiative  cooling  time  at  the  base.  Note  that  our 
b 

definition  for  X  has  a  sign  change  from  that  used  by  Field  (1965)  and  that 
the  negative  sign  is  required  for  stability.  Assuming  that  the 
free-surface  growth  rate,  ,  is  of  order  the  coronal  cooling  rate  1/t, 
Equations  (27)  and  (28)  yield  the  following  asymptotic  form  for  £* 
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£  ~  L 


,7/4 


-1/2 


[T0(0) 


where 


(1  +  T/Tb) 


(29) 
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The  boundary  conditions  must  now  be  specified.  Assuming  the  "natural" 
ones,  nlb  =  Tlb  =  0  leads  to: 


n  =  0  and  yfe  =  0 

at  each  end  of  the  loop.  Note  that  we  have  used  the  fact 
equilibrium  heat  flux  vanishes  at  the  base  in  order 
Equation  (23).  If  one  assumes  instead  the  rigid-wall 

V1  =  T  =  0,  this  leads  to  the  constraint  (Antiochos  1979): 

lb  lb 
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dx  =  0 


(23) 
that  the 
to  derive 
conditions 


(24a) 


and  the  conditions: 


+  yh5/u  n 

b  b 


=  0 


(24b) 


With  this  formulation  it  is  difficult  to  compare  the  effects  of  the  two 
sets  of  boundary  conditions  since  one  set  leads  to  homogeneous  conditions 
on  y,  (23),  whereas  the  other  does  not.  In  order  to  facilitate  the 

comparison,  let  us  for  the  moment  reformulate  the  problem  in  terms  of  a  new 

5/14  ,  .  „ 

variable  defined  as  4  =  y  +  y  n,  so  that  the  boundary  conditions  (24b) 

will  become  homogeneous.  Hence,  the  problem  reduces  to  solving  an  equation 
for  £  and  n  of  the  form  (14)  (only  the  form  of  G(x)  changes)  and  with 
boundary  conditions: 

for  free  surface,  t,  =0,  and  n  =  0; 

b  -1 

for  rigid  wall,  ?  =  0,  and  r,  =  /  4  y_5/14  dx  .  (25) 

•'o 

The  relation  between  the  two  sets  of  conditions  is  now  transparent.  We 
note  that  if,  as  expected,  c,  decreases  sufficiently  rapidly  in  the 
chromosphere,  then  (24)  implies  that  in  the  limit  of  large  chromospheric 
depths,  1. e. ,  L^/L  —  5  0, 

n  (Lc/L)  <C  y"5/14>  ,  (26) 


Under  the  special  assumptions  above,  the  perturbation  equations  (2)  - 
(6)  can  be  simplified  to  a  form  essentially  equivalent  to  the  one  we  used 
previously  (Antiochos  1979): 


dx  (y  dxU)  +  Q(x)VJ  +  X^  =  U  G(x)  , 


(14) 


where  we  have  defined  dimensionless  variables: 


=  s/L 


(15) 


L  is  the  loop  half-length, 


y(x)  =  (Tq/Tq (0) ) 


7/2 


(16) 


5/14  . 

=  -X _  y 

v  ds  1 


and  n  =  =  constant. 

0 


(17) 

(18) 


The  eigenvalue  X  is  defined  as: 


X  =  -  (L/Lcr  V  T  , 


(19) 


where  t  is  given  by  the  coronal  cooling  time  scale 


I  ’  >"6  Lc2  (Tq  (0)) 


-7/2 


(20) 


The  functions  Q(x)  and  G(x)  are  given  by: 

QM  =  -lilz  -  l  U  *2)2  +  aL  ZL) 

q  14  rfv2  y  \14  d*/  “  'T°  3T0' 


G(x)  =-y 


dx 
5/14 


3/5  A  +  dx^  +aT0  3tq  +aP0  8P0i 


(21a) 


(21b) 


where  a  is  a  constant  of  the  dimensions: 


a  =  106  L2  (T0(0))'772 


(22) 


T- - -  ^  I  m  m  V 


wall  where  the  pressure  P  jumps  back  up  to  its  coronal  value.  The 

lb 

physical  reason  for  this  is  straightforward.  Flows  in  the  corona  will 
cause  material  either  to  pile  up  or  to  evacuate  from  the  wall.  Hence, 
there  must  exist  a  pressure  gradient  near  the  wall  of  the  same  order  as  the 
coronal  gradient  in  order  to  decelerate  or  accelerate  the  flows.  However, 
over  most  of  the  chromosphere  the  flows  are  in  a  steady  state,  and  the 
pressure  perturbation  is  negligible. 

It  is  clear  from  these  results  that  the  appropriate  boundary  condition 
to  assume  at  the  loop  base  is  that  Plb  =  0.  This  condition  results 
naturally  in  the  chromosphere,  irrespective  of  what  i3  assumed  at  the  base 
as  long  as  the  base  is  placed  sufficiently  deep.  If  the  base  is  placed 
within  several  scale  heights  of  the  top  of  the  chromosphere,  then  it  is  the 
only  boundary  condition  that  is  appropriate;  the  use  of  rigid-wall 
conditions,  or  any  others,  is  likely  to  lead  to  incorrect  conclusions.  If 
the  base  is  placed  many  scale  heights  deep,  then  the  particular  boundary 
condition  used  becomes  irrelevant. 

III.  Analytic  Results 
a)  Perturbation  .Equations 

In  this  paper  we  calculate  growth  rates  and  eigenfunctions  only  for 
the  simplified  problem  discussed  by  previous  authors;  we  neglect  gravity 
and  any  spatial  variation  in  the  loop  area  or  heating  function.  We  discuss 
the  effects  of  gravity  and  the  area  variation  in  a  subsequent  paper 
(Antiochos  and  An  1985).  For  purposes  of  comparison,  we  calculate  the 
growth  rates  and  modes  for  several  sets  of  boundary  conditions,  although  as 
we  discussed  above,  the  physically  appropriate  ones  are  T^b  =  0  and 

"lb  = 


0 


P^Cs)  does  not 


perturbation  there  is  of  order  the  coronal  one  (i.e., 
diverge  exponentially  in  the  chromosphere).  For  a  base  that  is  several 
scale  heights  deep,  P  (0)/P_  is  negligible:  therefore,  the  form  of  the 

0  Ub 

coronal  perturbation  is  insensitive  to  the  boundary  condition  on  I’  or, 

equivalently,  V. 

lb 

The  situation  is  clarified  by  Figure  1,  where  for  a  given  T^s)  we 

plot  the  resulting  form  of  P^s)  and  r^V^s)  for  both  rigid-wall  Vlb  =  0 

and  free-surface  P  =0  boundary  conditions.  Since  this  is  for 

lb 

illustrative  purposes  only,  we  have  taken  an  extremely  simple  equilibrium 
loop  model  in  which  we  neglect  the  effects  of  gravity  in  the  corona  and 
assume  that  the  transition  region  is  vanishingly  thin.  The  corona  is, 

therefore,  isothermal  and  isobaric  and  is  defined  to  extend  from  s  =  0  to 

s  =  L_ .  Below  this  corona  is  an  isothermal  chromosphere  in  hydrostatic 
equilibrium  and  with  a  temperature  100  times  lower.  The  base  of  the  loop 
model  is  taken  to  occur  10  chromospheric  scale  heights  deep, 

L  =  Lc  +  10  Hb.  For  simplicity,  the  temperature  perturbation,  T^Tq,  is 
assumed  to  be  constant  in  the  corona  and  to  decay  exponentially  in  the 

chromosphere  on  a  size  scale  of  H  /10. 

b 

The  forms  for  P^(s)  and  nQVb(s)  obtained  from  (9)  and  (11)  are  plotted 
in  Figure  la  and  1b,  respectively.  We  have  taken  the  velocity  to  vanish  at 
the  top,  1^(0)  =  0,  which  would  be  the  case  for  the  lowest  order  mode  of  a 
symmetric  loop  model.  On  the  scale  of  this  figure,  the  modes  with  P  =  0 

and  those  with  V,,  =  0  are  indistinguishable  except  within  a  few  scale 
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heights  of  the  base.  In  this  boundary  region  the  rigid  wall  mass  flux 
exhibits  a  sudden  decrease  to  zero,  and  P.^  exhibits  a  sudden  increase  back 
up  to  its  coronal  value,  i.e.,  Plb  =  P^ ( 0 ) .  Note  that  the  magnitude  of 
this  jump  does  not  change  with  position  of  the  base.  No  matter  how  deep 


the  rigid  wall  is  placed,  there  will  occur  a  small  boundary  layer  at  the 


inodes  have  a  vanishingly  small  amplitude  in  this  region  compared  to  coronal 
amplitudes. 

(ii)  Rigid  Wall 

We  have  solved  Equation  (14)  for  the  rigid-wall  conditions  (24)  and 
for  the  exact  wide  range  of  equilibrium  models  used  with  the  free-surface 
conditions.  Some  of  the  results  are  shown  in  Figure  3.  The  growth  rate  of 
the  lowest  rigid-wall  mode  is  plotted  for  comparison  with  the  corresponding 
free-surface  mode.  The  rigid-wall  rate  is  seen  to  approach  the 
free-surface  one  as  the  chromospheric  depth  becomes  large,  1^,  /L  <<  1 ;  the 
ratio  of  the  growth  rates  is  near  unity,  vrA>j-  =  .85  for  the  deepest  model 
in  the  figure,  Lc/L  =  .07.  These  numerical  results  confirm  the  analytic 
arguments  above  that  for  large  depths  the  rigid-wall  condition  is 
equivalent  to  the  free-surface  one. 

We  find  that  the  rigid-wall  growth  rates  are  also  insensitive  to  the 

coronal  temperature,  the  base  temperature  or  the  form  of  the  heating. 

However,  Figure  3  shows  that  they  do  depend  on  the  depth  of  the  base 

region.  For  models  with  Lc  >  L/2  ,  the  growth  rate  is  negative,  indicating 

stability.  Note  that  this  value  for  Lc/L  actually  implies  a  large 

chromospheric  region.  Since  the  base  density  is  approximately  two  orders 

of  magnitude  greater  than  the  coronal  density,  the  bulk  of  the  loop  plasma 

is  in  the  base  for  Lc/L  =  .5.  Note  also  that  the  damping  rate  is  almost 
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constant  for  models  with  6  2  10  .  This  result  emphasizes  the  unphysical 

nature  of  the  rigid-wall  conditions  and  the  danger  in  using  them.  If  one 
assumes  the  rigid-wall  condition  and  then  considers  models  with  increasing 
chromospheric  depths,  it  first  appears  that  static  loops  are  stable  and 
that  the  stability  is  insensitive  to  base  conditions.  One  would  naturally 
conclude  that  the  models  are  physically  stable.  But  as  the  base  depth  is 


increased,  the  models  suddenly  become  unstable.  Somehow  the  addition  of 
too  much  stable  material  at  the  base  destabilizes  the  corona!  This  is 
clearly  an  unphysical  result  and  is  due  to  using  improper  assumptions  at 
the  onset. 

The  same  type  of  behavior  is  found  when  boundary  conditions  other  than 

T  =0  and  n  =0  are  used.  For  example,  Craig  and  McClymont  (1981)  find 

lb  lb 

that  the  conditions  dT  /ds  =  0  and  P.,  =  0  lead  to  abrupt  changes  in  the 

lb  lb 

growth  rate  at  very  small  chromospheric  depths.  We  have  also  examined 
these  conditions  and  a  variety  of  others  defined  by  the  vanishing  of  a 

linear  combination  of  T  and  F,,  or  P,,  and  V  .  Depending  on  the 

lb  lb  lb  lb 

particular  equilibrium  model,  the  boundary  conditions  assumed  can  result  in 
growth  rates  either  larger  or  smaller  than  the  free  surface  one;  however, 
in  all  cases  the  rates  tended  to  the  free-surface  result  for  large 
chromospheric  depths.  In  all  but  the  free-surface  case  the  growth  rates 
exhibited  an  abrupt  change  at  some  range  of  chromospheric  depth;  hence,  the 

conditions  T  =  n  =0  are  the  only  ones  that  are  generally  valid. 

lb  lb 

c)  Normal  Modes 

Along  with  the  eigenvalues  X  we  have  also  calculated  the 
eigenfunctions  y  for  all  the  equilibrium  models  considered.  TWo 
representative  cases  are  shown  in  Figures  4a  and  4b.  In  Figure  4a  we  show 
both  the  free-surface  Tj/Tq  and  the  rigid-wall  mode,  T-|A)T0»  for  the 

equilibrium  model  of  Figure  3  with  the  deepest  chromosphere  (L  /L  =  .07 

c 

corresponding  to  6  =  5  x  10  ~6).  In  Figure  4b  the  results  for  the  model 
with  6  =  10'3,  (Lc/L  =  .97),  are  shewn. 

We  first  discuss  the  results  for  the  deepest  model,  Fig.  4a.  In  order 
to  resolve  their  structure  in  the  transition  region,  the  modes  are  plotted 
against  a  nonlinear  distance  scale.  From  the  loop  apex,  s  =  0,  to  a  point 


two  coronal  depths  down,  s  =  2LC,  the  equilibrium  temperature  TQ  is  used  as 
the  independent  variable.  From  this  depth  down  to  the  loop  base,  s  =  L, 

the  distance  itself  is  used.  Of  course,  all  the  significant  structure  in 

the  eigenmodes  occurs  well  above  the  point  s  =  21.^.  The  region  belcw  this 

depth  is  included  in  the  figure  only  to  point  out  that  the  numerical 

calculations  verify  the  analytic  asymptotic  forms  obtained  above.  Using 
(29)  and  taking  x/xb  =0,  we  find  that  the  decay  scale  £  for  the 
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free-surface  mode  is  equal  to  2.16  x  10  L^.  Hence,  the  free-surface  mode 

should  have  essentially  zero  amplitude  in  the  deep  region.  The  numerical 
results  agree  with  this.  We  find  that  the  amplitude  of  T  /TQ  drops  belcw 

computer  accuracy  well  before  s  =  2L  .  On  the  other  hand,  the  amplitude  of 

the  rigid-wall  mode,  T^  /nTQ,  does  not  vanish.  Instead,  it  has  a  constant 
value  of  0.4  over  most  of  the  chromosphere  and  vanishes  only  near  the  base 

on  the  size  scale  £.  This  behavior  is  also  in  agreement  with  the  analytic 

results.  In  the  chromosphere,  T.^  /r\  TQ  =  %/n  ;  therefore,  from  (36) 

we  expect  the  mode  to  approach  the  value  2/5,  except  very  near  the  base, 
where  it  must  vanish  due  to  the  boundary  conditions. 

It  can  be  seen  from  Fig.  4a  that  the  modes  have  significant  amplitude 
only  in  the  corona  and  transition  region.  The  free-surface  mode  has  been 
normalized  to  equal  the  rigid-wall  value  at  the  apex.  Note  that  the 
rigid-wall  quantity  T^  /r\  Tq  does  not  have  an  aribtrary  normalization 
factor  since  this  quantity  is  the  ratio  of  the  temperature- to- pressure 
perturbations.  Its  magnitude  has  physical  significance.  Our  results 
indicate  that  the  pressure  perturbation  is  in  the  same  sense  but  smaller 
than  the  temperature  perturbation.  As  the  depth  of  the  chromospheric 
region  increases,  Lc/L  ->  0,  we  find  that  the  coronal  magnitude  of  T^/tiTq 
increases  as  (Lc/L)  \  This  confirms  our  previous  result  that  n  should 
decrease  as  Lc/L  for  a  fixed  temperature  perturbation. 


On  the  scale  of  Fig.  4a  the  two  modes  are  almost  identical  in  the 
corona  and  transition  region.  Only  near  the  base  is  there  a  significant 
deviation;  and  as  the  depth  of  the  base  region  is  increased,  the  two  modes 
become  indistinguishable.  The  modes  exhibit  a  stong  peak  near  5  x  It/1  K. 
At  approximately  this  temperature  the  radiative  losses  change  from  being  a 
destabilizing  to  a  stabilizing  effect,  so  that  belcw  this  temperature  the 
modes  decay  very  rapidly  on  the  scale  of  £.  It  is  important  to  emphasize 
that  although  the  modes  peak  strongly  in  the  transition  region,  the 
instability  is  not  a  local  one.  The  modes  have  finite  amplitude  throughout 
the  corona,  and  the  growth  rates  correspond  to  the  coronal  cooling  time 
rather  than  the  much  faster  cooling  time  in  the  transition  region.  Hence, 
the  instability  is  a  global  one  involving  the  complete  loop  except  for  the 
chromosphere. 

The  situation  is  somewhat  different  for  the  model  with  a  shallower 
base.  Fig.  4b.  The  free-surface  mode  is  very  similar  to  that  of  the  model, 
with  the  deeper  base.  Indeed,  this  mode  is  insensitive  to  just  about 
everything;  hence,  the  curve  shown  in  Fig.  4a  (or  4b)  can  be  considered  the 
universal  form. 

However,  the  rigid-wall  mode  is  obviously  very  different  from  that  of 

Fig.  4a.  The  most  striking  difference  is  that  the  mode  has  a  zero  crossing 

at  the  base  of  the  transition  region  and  is  negative  in  the  chromosphere. 

This  is  in  agreement  with  the  analytic  results.  We  note  from  (35)  that  as 

v  decreases  from  a  large  positive  quantity,  the  asymptotic  amplitude 

becomes  negative.  As  v  decreases  further  to  a  large  negative  quantity, 

the  amplitude  becomes  positive  again.  The  numerical  results  agree  with  the 
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behavior.  For  models  with  large  depths,  6  K  10  ,  the  amplitude  in  the 

chromosphere  is  positive.  For  models  with  intermediate  depths,  6  “10  , 

it  is  negative.  It  becomes  positive  again  for  small  depths,  6  >  10  .  For 


models  without  a  very  deep  base,  the  chromospheric  form  of  the  rigid-wall 
mode  is  clearly  quite  sensitive  to  the  exact  position  of  the  wall.  The 

free-surface  mode,  on  the  other  hand,  shews  no  change  for  all  values  of  6 

that  we  investigated. 

In  the  corona  the  rigid-wall  mode  has  the  same  general  shape  as  that 
in  Fig.  4a,  but  the  amplitude  is  an  order  of  magnitude  smaller.  This  means 
that  for  a  given  amplitude  temperature  perturbation,  the  effect  of  the 
decrease  in  the  base  depth  is  an  increase  in  the  coronal  pressure 

perturbation  by  an  order  of  magnitude.  In  fact,  the  pressure  perturbation 

actually  has  a  larger  coronal  amplitude  than  the  temperature  perturbation. 
We  note  from  Fig.  4b  that  the  ratio  of  the  temperature  to  pressure 
amplitude  is  less  than  unity  for  temperatures  down  to  ”  8  x  105  K.  This 
region  encompasses  about  70  percent  of  the  coronal  portion  of  the  loop. 
Since  from  (5) 


V 


1 


(44) 


our  results  imply  that  over  most  of  the  croona  the  temperature  and  density 
perturbations  must  be  in  the  same  sense.  A  decrease  (an  increase)  in  the 
temperature  is  accompanied  by  a  decrease  (an  increase)in  the  density  as 
well.  Since  a  decrease  in  the  temperature  implies  an  increase  in  the 
radiative  loss  rate  (for  T  >  105  K),  whereas  a  decrease  in  the  density 
implies  a  decrease  in  the  losses,  the  perturbations  oppose  each  other  as 
far  as  instability  is  concerned.  The  net  result  is  stability,  v  <0  for 
the  model  in  Fig.  4b.  For  models  with  a  deeper  base  region,  as  in  Fig.  4a, 
the  temperature  perturbation  is  much  larger  than  the  pressure  perturbation, 
and,  consequently,  the  density  perturbation  takes  the  opposite  sense.  Both 
the  temperature  and  density  disturbances  act  to  promote  instability,  and 
the  model  of  Fig.  4a  is,  indeed,  unstable,  v  >  0.  We  conclude  that  the 


stability  of  the  rigid-wall  models  such  as  Fig.  4b  is  a  boundary-condition 
effect.  It  is  due  solely  to  the  assumption  of  a  rigid-wall  at  a  parti eul a r 
depth. 

V.  Discussion 


There  are 

two  key 

results  of  this 

paper. 

The 

first  is 

that  the 

appropriate 

boundary 

conditions  to 

assume 

at 

the  loop 

base  are 

Tib  =  nlb  =  0.  For  models  with  a  sufficiently  deep  base 
(Lfe  »  min  (H^,!^)),  all  boundary  conditions  can  be  used;  however,  the  set 
above  is  the  only  one  that  is  appropriate  for  models  that  do  not  satisfy 
this  condition.  Indeed,  we  find  it  surprising  just  how  universally  valid 
these  boundary  conditions  appear  to  be.  Even  for  models  with  little  or  jqq 
chromosphere,  they  yield  the  correct  growth  rate  and  form  of  the  eigenmode 
for  the  instability.  Contrary  to  the  conclusions  of  several  authors 
(Habbal  and  Rosner  1979;  Chiuderi,  Einaudi,  and  Torricelli-Ciamponi 
1981;  Craig  and  McClymont  1981;  and  McClymont  and  Craig  198la,b,c),  the 
mere  presence  of  radiatively  stable  material  at  the  loop  base  does  not 
stabilize  the  static  models.  If  incorrect  boundary  conditions  are  used, 
then  the  base  can  act  either  to  stabilize  or,  as  we  showed  above,  to 
destabilize  the  static  models.  However,  if  the  proper  conditions  are  used, 
then  the  amount  of  base  material  has  no  effect  on  the  stability. 

This  is  a  convenient  result  from  a  numerical  point  of  view.  It  is 
often  quite  difficult  to  handle  numerical  models  with  extremely  deep  base 
regions  (e.g. ,  Craig,  Robb,  and  Rollo  1982).  Our  results  indicate  that 
models  with  very  shallow  bases  can  yield  correct  results  if  the  correct 
boundary  conditions  are  used. 


The  second  main  result  of  this  paper  is  that  under  our  assumptions, 


static  models  of  coronal  loops  are  thermally  unstable.  The  growth  rate  of 
the  instablity  is  of  the  order  the  coronal  cooling  time  (20).  The  basic 
assumptions  that  we  made  for  the  equlibrium  models  were  that  gravity  and 
all  spatial  variations  in  the  loop  area  or  heating  rate  are  negligible. 
These  are  the  most  commonly  used  assumptions  in  calculating  static  models. 
An  important  question  is:  hew  sensitive  are  our  results  to  these 
assumptions?  McClymont  and  Craig  (1981a,  b,c)  have  considered  the  question 
of  spatial  variations  in  the  heating  rate.  They  find  that  these  can  act 

either  to  enhance  or  to  damp  instability.  Lacking  a  model  for  the  coronal 

heating  process,  nothing  definitive  can  be  stated  on  this  issue;  hence,  we 
will  not  consider  it  further.  We  will  discuss  the  effects  of  gravity  and 
the  area  variation  in  a  forthcoming  appear.  In  general,  the  area  variation 
acts  to  enhance  the  instability;  however,  we  do  not  expect  it  to  enhance 
the  growth  rates  to  values  much  beyond  the  ones  here,  since  these  are 
already  of  the  order  of  the  coronal  cooling  rate.  Gravity,  on  the  other 
hand,  acts  to  damp  the  instability  (e.g. ,  Wragg  and  Priest  1982);  hence, 
our  results  are  certain  to  apply  only  for  loops  in  which  gravity  can  be 
neglected. 

We  can  identify  at  least  two  possible  physical  situations  in  the  solar 
corona  where  the  effects  of  gravity  should  be  small  so  that  our  results 
clearly  apply.  One  is  the  case  of  loops  sufficiently  lew  lying,  heights 
K  1000  km,  so  that  the  loop  height  is  smaller  than  any  gravitational  scale 
height  in  the  loop.  Another  is  the  case  of  magnetic  field  lines  that  are 
concave  in  the  corona,  so  that  near  the  loop  center  the  force  of  gravity  is 
directed  toward  the  center  rather  than  away  from  it,  as  in  the  usual  case. 
This  is  the  type  of  magnetic  geometry  believed  to  be  responsible  for 

quiescent  prominences.  Our  results  imply  that  these  two  types  of 

structures  are  naturally  unstable  without  any  special  requirements  on  the 


coronal  heating  process.  We  believe  that  the  thermal  instability  discussed 
in  this  paper  is  responsible  for  the  cool  material  that  is  often  observed 


to  occur  in  the  corona. 

Of  course,  in  order  to  prove  this  statement,  the  nonlinear  development 
of  instability  must  be  calculated.  Clearly,  the  instability  must  not 
saturate  at  a  lew  amplitude  if  it  is  to  produce  observational  effects  in 
the  corona.  Only  the  linear  growth  rates  were  considered  in  this 
paper;  the  nonlinear  evaluation  must  be  determined  by  numerical  simulation. 
Numerical  simulations  of  coronal  loops  have  been  performed  by  several 
authors  (e.g. ,  Craig,  Robb,  and  Rollo  1982;  Peres  ££  al.  1982;  Oran, 
Mariska,  and  Boris  1982).  These  simulations  have  not  exhibited  any 
significant  evidence  for  thermal  instability;  however,  to  our  knowledge  no 
one  has  considered  models  in  which  the  effects  of  gravity  are  negligible 
and  the  correct  boundary  conditions  are  employed  at  the  loop  base.  Our 
results  indicate  that  the  instability  should  be  present  in  such  models.  We 
intend  to  investigate  the  nonlinear  stability  of  coronal  loops  in  a  future 
work. 
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Figure  Captions 


Figure  1  Plots  of  free-surface  and  rigid-wall  pressure  perturbation 
(Figure  la)  and  mass  flux  perturbation  (Figure  1b)  for  the 
highly  simplified  model  described  in  the  text.  The 
free-surface  results  are  indicated  by  solid  curves  and  the 
rigid-wall  one  by  broken  curves;  the  two  types  of  modes  are 
distinguishable  only  near  the  base,  s  =  L. 

Figure  2  The  radiative  loss  coefficient  of  Raymond,  Cox  and  Smith  (1976) 

(broken  line)  and  our  analytic  fit  to  it  (solid  line).  The  units 
in  the  figure  are  arbitrary. 

Figure  3  The  free-surface  (solid  curve)  and  rigid-wall  (broken  curve) 
growth  rates  for  a  series  of  static  models  with  varying  base 
temperatures  and/or  varying  base  depths. 

Figure  4a  Plot  of  lcwest  free-surface  temperature  perturbation  (T^/Tq) 

and  rigid-wall  perturbation  (T^Pq/TqP^)  for  the  e9uilibrium  model 
with  6  =  5  x  10-6  and  a  length  L  =  14.2  Lc<  The  free-surface 
mode  is  indicated  by  a  solid  curve,  and  the  rigid-wall  one  by  a 
broken  curve;  they  are  indistinguishable  in  the  corona  and 
transition  region.  The  scale  of  the  abscissa  is  log  (TQ)  for 

0  i  s  i  2L  ,  and  s/L  for  2L  <  s  i  L. 

c  c 

-3 

Figure  4b  Same  as  in  Figure  4a  but  for  the  equlibrium  model  with  6=10  , 

and  with  an  abscissa  of  log  (TQ)  for  0  £  s/L  <(  .94  and 
log  (s/L)  tor  .94  i  s/L  <  1.  Note  that  the  rigid-wall 
amplitude  is  negative  in  the  region  s/L  >  .94. 
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